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ABSTRACT 

The neutrino mechanism of core-collapse supernova is investigated via non-relativistic, two- 
dimensional (2D), neutrino radiation-hydrodynamic simulations. For the transport of electron flavor 
neutrinos, we use the interaction rates defined by Bruenn (1985) and the isotropic diffusion source 
approximation (IDSA) scheme, which decomposes the transported particles into trapped particle and 
streaming particle components. Heavy neutrinos are described by a leakage scheme. Unlike the 
“ray-by-ray” approach in some other multi-dimensional supernova models, we use cylindrical coordi¬ 
nates and solve the trapped particle component in multiple dimensions, improving the proto-neutron 
star resolution and the neutrino transport in angular and temporal directions. We provide an IDSA 
verification by performing ID and 2D simulations with 15 and 20 M 0 progenitors from Woosley et 
al. (2007) and discuss the difference of our IDSA results with those existing in the literature. Ad¬ 
ditionally, we perform Newtonian ID and 2D simulations from prebounce core collapse to several 
hundred milliseconds postbounce with 11, 15, 21, and 27 M 0 progenitors from Woosley et al. (2002) 
with the HS(DD2) equation of state. General relativistic effects are neglected. We obtain robust 
explosions with diagnostic energies i^dig > 0.1 — 0.5 B for all considered 2D models within approx¬ 
imately 100 — 300 milliseconds after bounce and find that explosions are mostly dominated by the 
neutrino-driven convection, although standing accretion shock instabilities are observed as well. We 
also find that the level of electron deleptonization during collapse dramatically affect the postbounce 
evolution, e.g. the ignorance of neutrino-electron scattering during collapse will lead to a stronger 
explosion. 

Subject headings: hydrodynamics — instabilities —neutrinos — supernovae: general 


1. INTRODUCTION 

After nearly five decades from the first attempt to 
obtain a neutrino-driven explosion by Colgate & White 
(1966), the explosion mechanism of core-collapse super¬ 
novae (CCSNe) remains elusive. It is generally believed 
that neutrino transport and convection are important in¬ 
gredients to achieve a successful explosion (see recent re¬ 
views in Janka 2012; Burrows 2013; Foglizzo et al. 2015). 
However, modeling CCSNe with Boltzmann transport 
in three dimensions is numerically expensive and too 
time consuming with current computing resources, since 
the neutrino radiation hydrodynamics with Boltzmann 
transport is essentially a seven-dimensional problem: 
three spatial coordinates, two angular degrees of free¬ 
dom, neutrino energy and time. In addition, there are 
three types of neutrino species and their anti-particles 
that would require a solution of the Boltzmann equation. 

CCSN simulations with Boltzmann transport have 
been studied in ID (Mezzacappa & Bruenn 1993; 
Liebendorfer et al. 2001, 2004; Sumiyoshi et al. 2005), 
in 2D (Livne et al. 2004; Ott et al. 2008; Brandt et al. 
2011) and recently in 3D with low resolutions and 
fixed background profiles (Sumiyoshi et al. 2015). How¬ 
ever, these 2D and 3D works have ignored the velocity- 
dependent terms and decoupled these from the energy 
groups for simplicity, and the spatial resolutions are not 
sufficient to achieve a correct turbulence cascade. There¬ 
fore, approximate methods for the neutrino transport in 
multi-dimensional simulations are still necessary at this 
moment. 


Simple approximatation schemes include the light- 
bulb scheme (Murphy & Burrows 2008; Hanke et al. 
2012; Couch 2013), neutrino leakage (Janka & Mueller 
1996; Rosswog & Liebendorfer 2003; O’Connor & Ott 
2011; Couch & O’Connor 2014), and gray transport 
(Fryer et al. 1999; Scheck et al. 2006). In these schemes, 
the neutrino transport is relatively cheap and therefore, 
it is possible to perform high resolution simulations in 
3D (with effective angular resolutions < 1°), better de¬ 
scribing the turbulence and convection behind the shock 
and the standing accretion shock instabilities (SASI, 
Blondin et al. 2003). However the neutrino transport in 
these schemes is possibly over-simplified because it does 
not really follow the transport of the neutrino distribu¬ 
tions and therefore, cannot describe the neutrino heating 
self-consistently. 

A more sophisticated but still approximated scheme 
for the neutrino transport is called the Moment scheme: 
The multi-group flux limited diffusion (MGFLD) scheme 
(Bruenn et al. 2013; Dolence et al. 2015) takes the ze¬ 
roth angular moment of neutrino moment expansions 
(i.e. energy density). The Ml moment scheme addi¬ 
tionally evolves the momentum density and considers 
the higher moment closure with an analytic form 
(Pons et al. 2000; Kuroda et al. 2012; O’Connor & Ott 
2013; Obergaulinger et al. 2014; O’Connor 2014; 
Kuroda et al. 2015) or by a variable Eddington 
tensor (Burrows et al. 2000; Rampp & Janka 2002; 
Thompson et al. 2003; Buras et al. 2006; Muller et al. 
2012; Melson et al. 2015). However, as reported by 
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Kuroda et al. (2015), the Ml scheme in 3D general 
relativity (GR) is still very expensive and difficult to 
run in long term postbounce simulations with high 
resolutions (an effective angular resolution < 4°). 

Another approximated transport scheme is the 
isotropic diffusion source approximation (IDSA, 
Liebendorfer et al. 2009). In the IDSA, the distribution 
function is separated into an opaque trapped particle 
component and a transparent streaming particle com¬ 
ponent. The two components are linked by a source 
term. Therefore, the transport equation becomes a 
diffusion problem in the opaque region and a ray tracing 
problem in the transparent region (Liebendorfer et al. 
2009). Multi-dimensional simulations with the IDSA 
have been performed by Suwa et al. (2011, 2014) in 
2D, and Takiwaki et al. (2014, 2012) in 3D. However, 
most of these multi-dimensional studies with the IDSA 
or with other transport schemes do not really solve 
the neutrino transport in multi-dimensions. Instead, 
they consider the angular distribution to be several ID 
problems, i.e. apply the ray-by-ray (RbR) approach. 
Sumiyoshi et al. (2015) and Dolence et al. (2015) have 
pointed out that the RbR approach may artificially 
exaggerate the neutrino flux variations in the angular 
and temporal components, since the temporal variation 
of the neutrino fluxes in convective regions is greatly 
ignored in the RbR approach. Additionally, the RbR 
approximation overestimates the heating of accretion 
luminosity on its own ray and underestimates the 
heating in neighbor rays. 

In this Paper, we present two-dimensional CCSN sim¬ 
ulations with the IDSA for neutrino transport in cylin¬ 
drical coordinates, which is in principle easy to extend to 
3D and has better resolution and boundary conditions for 
the description of the proto-neutron star (PNS) than sim¬ 
ulations in spherical coordinates. We solve the trapped 
particle component in multi-dimensions to improve the 
neutrino transport in angular and temporal directions. 
Furthermore, we study the neutrino transport effects 
during collapse by comparing our IDSA solver with a 
parametrized deleptonization scheme from Liebendorfer 
(2005). We find that the postbounce explosion dynamics 
is sensitive to the detailed neutrino interactions before 
core bounce, such as neutrino-electron scattering (NES) 
and electron-capture rates. Additional effects such as 
GR and magnetic fields are also crucial factors on study¬ 
ing CCSN explosion mechanism, in particular for fast- 
rotating and more massive progenitors, but we leave 
these parts for future work. In the following section, the 
numerical methods and the IDSA implementations are 
described. A verification of our IDSA implementation 
and a comparison with other neutrino transport schemes 
is shown in Section 3. In Section 4, we apply our IDSA 
implementation to multiple progenitors, using a new SN 
equation of state (EOS). In Section 5, we summarize our 
results and conclude. A discussion of the different EOS 
is presented in Appendix A. 

2. NUMERICAL METHODS AND MODELS 

We describe the numerical code and the correspond¬ 
ing setup of our simulations in 2.1. The method and 
implementation of the IDSA for neutrino transport is 
demonstrated in 2.2. We present the investigated EOS 
and supernova progenitors in Sections 2.3 and 2.4. 


2.1. Numerical Code And Initial Setup 

We use FLASH 1 version 4 (Fryxell et al. 2000; 
Dubey et al. 2008; Lee 2013) to solve the Eulerian equa¬ 
tions of hydrodynamics, 



I + v-M-o 

(i) 

- + V • (pvv) + VP = -pV<f> 

( 2 ) 

dpE 

dt 

+ V • [(pE + P)v\ = —pv • 

( 3 ) 


^El+V-(pvY e )=0 

( 4 ) 


^ + v-K)= o 

( 5 ) 

7 + V • ( v(pZf A 4 ) = 0 

( 6 ) 


where p is the gas density, v is the velocity vector, P is 
the gas pressure, E is the total specific energy, <f> is the 
gravitational potential, Y e is the electron fraction, Y l is 
the particle number fraction, and Z is the particle mean 
specific energy. The index l labels different species of 
trapped particles (i.e. v e and z7 e ). A detailed description 
of the neutrino-transport method will be presented in 
Section 2.2. 

FLASH is a parallel, multidimensional hydrodynamics 
code based on block-structured Adaptive Mesh Refine¬ 
ment (AMR). Our simulation setup is essentially simi¬ 
lar to what has been implemented in Couch (2013) and 
Couch & O’Connor (2014), but replacing their radiation 
transfer solver by an IDSA solver (Liebendorfer et al. 
2009). We use the third-order piecewise parabolic 
method (PPM, Colella & Woodward 1984) for spatial re¬ 
construction, the HLLC Riemann solver and the Hybrid 
slope limiter for shock-capture. The “hybrid” Riemann 
solver is widely used to avoid an odd-even numerical in¬ 
stability when the shock is aligned with the grid (Quirk 
1991). However, we don’t see this instability in our 
CCSN simulations by comparing simulation results with 
the hybrid Riemann solver. The HLLC Riemann solver 
shows a better turbulent cascade based on the implicit 
large eddy simulations (Radice et al. 2015). Effects from 
general and special relativity and from magnetic fields 
are ignored. 

Simulations are performed in ID spherical and 2D 
cylindrical coordinates. The center of a progenitor star is 
located at the origin of the simulation box. The simula¬ 
tion box includes the inner 10 4 km in radius of a progeni¬ 
tor star in ID and the full 180° in cylindrical coordinates 
with r max = z max = 10 4 km and z m i n = —10 4 km in 2D. 
The central r < 100 km sphere has the smallest zone 
size of 0.488 km and the AMR mesh is dynamically ad¬ 
justed based on the magnitude of the second derivatives 
of gas density, pressure, and entropy. To save computa¬ 
tion time, we apply additional AMR decrements based on 
the distance to the origin. E.g. the first AMR decrement 
is enforced at r ~ 100 km and the second AMR decre¬ 
ment at r ^ 200 km, the next at r ~ 400 km, and so 

1 http://flash.uchicago.edu 
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Table 1 

Neutrino Interactions included in the Simulations 


Neutrino Interaction Reference 


1. v e +n ^ e~ +p Bruenn (1985) 

2. u e p ^ e + + n Bruenn (1985) 

3. v e + A ^ e~ + A' Bruenn (1985) 

4. v + a/A/N v + ot/A/N Bruenn (1985) 

5. v + e ± ^ v + Mezzacappa & Bruenn (1993); Liebendorfer (2005)’!' 


Note. — n= free neutrons, p= free protons, N= free neutrons or protons, A= nuclei 
besides a particles, a=alpha particles, u e = electron type neutrinos, r> e = anti-electron type 
neutrinos, u= all type of neutrinos, e _ =electron, and e“*~ = positron. 

^ In this Paper, we include the neutrino-electron scattering by means of using the 
parametrized deleptonization in Liebendorfer (2005) in the pre-bounce phase, since this 
reaction is only important during core collapse. 


on. The maximum zone size is 62.5 km. We employ the 
“outflow” boundary condition at the outer boundaries 
and the “reflect” boundary condition at the inner bound¬ 
aries. The gravitational potential is solved by the new 
improved multipole solver in FLASH (Couch et al. 2013) 
with a maximum Legendre order, / max = 16. 

It should be noted that the “outflow” boundary condi¬ 
tion, as reported by Couch (2013), causes a zero-gradient 
mass accretion at the boundary which will overestimate 
the inflow and suppress the explosion at late time. We 
therefore extend our simulation domain to 10,000 km to 
minimize the effect of boundary conditions. 

2.2. Neutrino Transport 

In the IDSA (Liebendorfer et al. 2009; Berninger et al. 
2013), we decompose the distribution function / of trans¬ 
ported particle species and the neutrino transport opera¬ 
tor D into a trapped particle component and a streaming 
particle component. We assume that these two com¬ 
ponents evolve separately. With this assumption, we 
rewrite the transport equation, D(f = f l + f s ) = (7, 
where C = C l + C s is a collision integral, by two equa¬ 
tions, 

D(f) = C* - E (7) 

D(f s ) =C S + E (8) 

where E is the diffusion source term, which converts 
trapped particles (t) into streaming particles (s) and vice 
versa. 

By using the diffusion limit, the diffusion source term 
E in the trapped transport equation could be expressed 
by (see Liebendorfer et al. 2009 for a more detailed dis¬ 
cussion), 

E = min | max 
where, 

« = v -(5(JT^) V 4 (10> 

is a non-local diffusion scalar, j is the emissivity, x is the 
absorptivity, </> is the opacity, and /a is the angle cosine. 
The distribution functions / t and f s can be solved by us¬ 
ing equations (8) and (10) in Liebendorfer et al. (2009). 
Once the distribution functions are known, the net in¬ 
teraction rates of transported particles can be calculated 


and further updates of the fluid quantities, such as fl, Y e , 
E, Yf, and Zj, can be derived. 

We solve the streaming component in ID with the orig¬ 
inal IDSA solver, but solve the trapped component and 
a in multi-dimensions. In each time step, we take the 
angle-averaged radial bins of fluid and neutrino quan¬ 
tities, such as p, s, T, y e , Yf, as ID inputs for 
the streaming component. The radial bins contain 600 
zones sampled from the center of the proto-neutron star 
(PNS), defined by the location of maximum density, up 
to r = 5,000 km. The radial bins are uniformly spaced 
(Ar = 2 km) up to a radius of r = 400 km. Be¬ 
yond 400 km, the zone spacing logarithmically increases. 
The streaming component carries only to the location 
of neutrino spheres R U (E ) and the angular integration 
of the spectral neutrino flux F S (E). The local heating 
from streaming neutrinos are then determined based on 
the local neutrino interaction rates based on the multi¬ 
dimensional hydro quantities. 

It should be noted that the assumption of spherical 
symmetry in the streaming component may calculate 
the neutrino flux and heating incorrectly when shocks 
become highly asymmetric, especially in 2D simulations 
after the SASI has developed. On the other hand, a 
RbR approach may artificially enhance the asymmetry 
and lead to incorrect results as well (Sumiyoshi et al. 
2015; Dolence et al. 2015). It is important to have multi¬ 
dimensional simulations with both approximations be¬ 
cause they are likely to bracket the actual solution: The 
RbR focusses all heating from the accretion luminosity 
on its own ray, while the angular integration distributes 
it over all directions. In reality, a solution in between 
these extremes is expected. 

For the trapped component, we explicitly solve the dif¬ 
fusion source E and a in multi-dimensions and update f l 
locally. Together with the streaming component, the new 
multi-dimensional interaction rates can be updated. The 
new neutrino sources are re-evaluated based on multi¬ 
dimensional neutrino sources and production rates, and 
then used for the next streaming step. Since we solve 
the diffusion part explicitly, a stable solution requires 
a small neutrino time step t v ^ Ax 2 /(2Xc) ~ Ax/2c. 
Therefore, in our ID and 2D simulations, we use a fixed 
hydro timestep with thydro = 10 -6 sec and do sub-cycling 
for the neutrino transport. Since the IDSA solver in the 
streaming component allows larger time step, we adopt 
two types of sub-cycles. One for the streaming compo- 


a - 


(j + J f s d^0 ,j\ 


( 9 ) 
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Table 2 

Parameters for the Parametrized Deleptonization 


Progenitor 

sll.O (W02) 

sl5.0 (W02) 

s21.0 (W02) 

s27.0 (W02) 

sl5 (W07) 

s20 (W07) 

Pi [g cm -3 ] 

1.5 x 10 8 

9.0 x 10 7 

2.5 x 10 8 

2.5 x 10 8 

2.2 x 10 8 

3.0 x 10 8 

p 2 [g cm- 3 ] 

1.2 x 10 13 

9.0 x 10 12 

1.0 x 10 13 

1.0 x 10 13 

9.5 x 10 12 

1.0 x 10 13 

Vi 

0.5 

0.5 

0.5 

0.5 

0.5 

0.5 

v 2 

0.287 

0.282 

0.279 

0.279 

0.279 

0.273 

v c 

0.02 

0.03 

0.017 

0.017 

0.022 

0.017 


Note. — Parameters used in the fitting formula from Liebendorfer (2005) 


nent with K = 5x 10 7 sec, and the other for the trapped 
component with t l v = 10 _r sec. Additionally, the neu¬ 
trino pressure contributes extra momentum and can be 
expressed by, 


dv = _l fpZl\ 
dt p \3rab ) 


(ii) 


We note that Equations (6) and (11) include all 
thermodynamically important 0(v/c ) terms of the 
Boltzmann transport equation. These 0{y/c ) terms 
have been considered as crucial for CCSN modeling 
(Mezzacappa & Bruenn 1993; Liebendorfer et al. 2009; 
Lentz et al. 2012a). 

Our IDSA solver only includes electron flavor neu¬ 
trinos. We use 20 energy bins which are spaced 
logarithmically from 3 MeV to 300 MeV. For heavy 
neutrinos, the energy release is treated by a leak¬ 
age scheme that is based on the local diffusion and 
the local production rates (Hannestad & Raffelt 1998; 
Rosswog & Liebendorfer 2003; Perego et al. 2014). 

Table 1 summarizes all weak interactions that are in¬ 
cluded in this work. All weak interactions in the IDSA 
solver are using the Bruenn description (Bruenn 1985). 
Note that the IDSA solver dose not include the NES 
(interaction 5 in Table 1). NES is important during 
the collapse phase but provides minor contributions in 
the postbounce phase (Bruenn 1989) (see Section 3.2 for 
a demonstration). Liebendorfer (2005) found a simple 
relation between the electron deleptonization with the 
gas density in the collapse phase, where electron fraction 
and entropy can be parametrized by density and chem¬ 
ical potentials. This parametrized deleptonization (PD) 
scheme can take interactions effectively into account that 
are only implemented in the Boltzmann solver that is 
used to determine the parameters for the PD. We have 
recalibrated the fitting parameters with different progen¬ 
itors for the PD scheme by running AGILE-BOLTZTRAN 
(Liebendorfer et al. 2004), since the original parameters 
in Liebendorfer (2005) are calibrated for the progeni¬ 
tors G15 (sl5s7b2, Woosley & Weaver 1995) and N13 
(Nomoto & Hashimoto 1988). Table 2 summarizes the 
fitting parameters we use during collapse. In this Paper 
we include the effect of NES by means of using the PD 
scheme in the collapse phase. Since the scheme is inde¬ 
pendent of the neutrino fractions, we update the neutrino 
fractions and energy through the IDSA solver during col¬ 
lapse. After bounce, we turn off the PD and switch to 
the IDSA solver. For simulations without NES, we use 
the IDSA solver from the beginning of core collapse to 
postbounce. 

To verify our IDSA implementation in FLASH, we pro¬ 


vide a code comparison of FLASH-IDSA with AGILE-IDSA 
and with existing results from the literature in Section 3. 

2.3. Equation of State 

We use the nuclear EOS unit in FLASH which in¬ 
corporates the finite temperature EOS routines from 
O’Connor & Ott (2010) and Couch (2013) 2 . The Lat- 
timer & Swesty EOS (with incompressibility, K = 
220 MeV; Lattimer & Swesty 1991) and Hempel & 
Schaffner-Bielich (HS) DD2 EOS are used in this work. 
The HS(DD2) EOS employs the density-dependent rel¬ 
ativistic mean-field interactions of Typel et al. (2010). 
The description of nuclei in supernova matter is based 
on (Hempel & Schaffner-Bielich 2010). This EoS was 
first applied in core-collapse supernova simulations by 
Fischer et al. (2014), where further details can be found. 
LS220 is one of the most common and well-studied 
EOS in the supernova community. However, it has 
some deficiencies, for example it is based on the sin¬ 
gle nucleus approximation for heavy nuclei, and con¬ 
siders only the alpha particle as a degree of freedom of 
all possible light nuclei. See Hempel et al. (2015) for a 
comparison of predictions for cluster formation for the 
HS(DD2) and the LS220 EOS with experiments, where 
good agreement was found. Furthermore, it was shown 
by Kruger et al. (2013); Fischer et al. (2014), that the 
neutron matter EOS of LS220 is in disagreement with 
constraints from Chiral effective field theory. Further¬ 
more, its low-density symmetry energy deviates from 
constraints obtained from finite nuclei, see Fig. 9 of 
Hempel (2014). No multi-dimensional simulations have 
been performed with HS(DD2) at this moment. We use 
LS220 for a code verification test in Section 3 and then 
use DD2 for our main simulations in Section 4. A low- 
density extension for the EOS is included in the routines 
from O’Connor & Ott (2010). In Appendex B, we pro¬ 
vide a brief discussion of the differences between LS220 
and DD2. 


2.4. Supernova Progenitors 

In this Paper, we consider four different nonrotat¬ 
ing, solar-metallicity progenitors, si 1.0, sl5.0, s21.0, 
and s27.0 from Woosley et al. (2002) 3 for our multi¬ 
ple progenitor study. We also consider two nonro¬ 
tating, solar-metallicity progenitors, sl5 and s20 from 
Woosley & Heger (2007) for a comparison study. Fig¬ 
ure 1 shows the initial density distribution of these six 
progenitors. The si 1.0 progenitor has the highest core 

2 http://www.stellarcollapse.org 

3 http://2sn.org/stellarevolution/ 
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Figure 1. Density as function of enclosed mass for the four 
investigated progenitor models from Woosley et al. (2002) and 
Woosley & Heger (2007). 

density but the most dilute envelope. The s21.0 and 
s27.0 and s20 progenitors have a similar density distri¬ 
bution and the most massive envelope among the mod¬ 
els, but have a different mass of the iron core and Si/O 
shell. The locations of regions with a high density gra¬ 
dient correspond to the Si/O interface. For the same 
progenitor mass, si5 has a denser core and more mas¬ 
sive envelope than si5.0. Another common progenitor 
model used in the literature is the sl5s7b2 progenitor 
from Woosley & Weaver (1995). It has the same pro¬ 
genitor mass as s 15.0 and si5 but the Si/O interface is 
located much further inside. This difference may make 
significant changes on the postbounce shock-radius evo¬ 
lutions when the shock reaches the interface due to a 
different mass-accretion history (Suwa et al. 2014). 

To adopt the progenitor models in FLASH, we map the 
one-dimensional density, temperature, electron fraction, 
and radial velocity profiles from Woosley et al. (2002) 
into our 1D/2D grids in FLASH. Other thermodynamic 
quantities are re-calculated using the EOS solver in 
FLASH. Neutrino fractions are set to zero at beginning. 

3. IDSA VERIFICATION 

To verify the IDSA implementation in FLASH, we 
first compare our ID FLASH simulations with simu¬ 
lations with AGILE-IDSA (Liebendorfer et al. 2009) in 
Section 3.1. AGILE-IDSA is a ID spherically sym¬ 
metric Lagrangian code which is publicly available on¬ 
line 4 . In Liebendorfer et al. (2009), a nice agree¬ 
ment of AGILE-IDSA with the GR Boltzmann code 
AGILE-BOLTZTRAN (Liebendorfer et al. 2004) was shown. 
Since we want to compare our results with the same neu¬ 
trino transport scheme (IDSA), we run additional sim¬ 
ulations in AGILE-IDSA but turn off the GR correction 

4 https://physik.unibas.ch/~liebend 


Figure 2. ID radial profiles in mass coordinate during the col¬ 
lapse of the sl5.0 progenitor with LS220 EOS. Different colors 
represent different times when the central density reaches p c = 
10 11 ,10 12 ,10 13 ,10 14 and at bounce. Solid lines show simulation 
results from FLASH and dashed lines show simulation results from 
AGILE-IDSA. Note that there are some slight time offsets between 
the results from FLASH and AGILE-IDSA because their output hies 
do not exactly match the given central densities. 

for the gravitational potential. We also turn off the PD 
during collapse in both codes to test the IDSA solver for 
neutrino transport before and after bounce. 

In Section 3.2, we show ID and 2D FLASH simulations 
with the two widely used progenitors, si5 and s20 from 
Woosley & Heger (2007) and discuss the differences com¬ 
pared to other transport schemes. Furthermore, we run 
AGILE-BOLTZTRAN simulations to demonstrate the impor¬ 
tance of NES during core collapse and discuss the influ¬ 
ence of the utilized neutrino opacities and electron cap¬ 
ture rates on the shock evolution and neutrino signal. 

3.1. Code Comparison with AGILE-IDSA 

We perform ID CCSN simulations of the four inves¬ 
tigated progenitors (sll.0, 515.0, 521.0, and 527.0) with 
the LS220 EOS in both FLASH and AGILE-IDSA. Figure 2 
shows the radial profiles of density, electron fraction, en¬ 
tropy and radial velocity of the progenitor 515.0 during 
collapse. It is clear to see that the bounce profiles (black 
lines in Figure 2) are nearly identical in the two codes. 
Small differences at the center could originate from the 
fact that AGILE-IDSA is a Lagrangian code with a mov¬ 
ing mesh. At the beginning of a simulation, the inner¬ 
most region of the progenitor star is much less resolved 
than with the Eulerian grid of FLASH. Therefore, the ra¬ 
dial profiles in AGILE-IDSA evolve slightly slower than 
those in FLASH. Furthermore, the bounce time in FLASH 
is about 0.1 — 12 ms later than in AGILE-IDSA, depending 
on the progenitor star. 

Figure 3 shows the electron and lepton fractions of the 
progenitor 515.0 at core bounce. It is found that the elec¬ 
tron fraction is consistent in both codes, but lepton frac- 
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Figure 3. Electron (solid lines) and lepton (dashed lines) fraction 
profiles at core bounce in model sl5.0. Red lines show simulation 
data from Agile-IDSA and blue lines represent data from FLASH. 
The black line indicates the electron fraction from the fitted PD 
scheme (Table 2). Both Agile-IDSA and FLASH are ID simulations 
with LS220 EOS and without NES. 

tion shows a mismatch at 10 10 < p < 10 12 g cm -3 , sug¬ 
gesting a lower electron neutrino distribution at low den¬ 
sities region in FLASH. However, this difference does not 
alter the postbounce evolution much since the neutrinos 
are simply free streaming in the low density region. We 
note that in our multi-dimensional IDS A solver, we have 
employed a limiter for the diffusion scalar a. The limiter 
enforces the trapped component to have a small value, 
when r > 1.25 x R U (E ), to prevent unphysical oscilla¬ 
tion on the <a, when the trapped neutrino density goes 
to zero in the free-streaming regime. This limiter leads 
to small differences of the neutrino luminosity and mean 
energy, but does not affect the hydrodynamic quantities. 
Figure 4 shows the particle spectra of the trapped and 
streaming components at 150 ms postbounce of the pro¬ 
genitor si5.0. We show three different regions where the 
trapped particle component dominates (at r = 40 km), 
trapped and streaming particle components are compa¬ 
rable (at r = 60 km), and where the streaming particle 
component dominates (r = 100 km). There are small 
differences in the particle spectra, but we do not expect 
exactly identical spectra in the two codes, since the hy¬ 
drodynamic parts are different. 

A comparison of radial profiles of the density, elec¬ 
tron fraction, entropy, and radial velocity of the progen¬ 
itor si5.0 in both codes at different postbounce times 
is shown in Figure 5. We note that although there are 
some differences in the radial shock location, the profiles 
and shock locations are consistent in the mass coordi¬ 
nate. The post bounce shock radius evolution for the 
four progenitor models are shown in Figure 6. Overall, 
the evolutions of the shock radius are nearly the same 
but show slight differences at ^ 120 ms (~ 220 ms) in 


Figure 4. Particle spectra at 150 ms postbounce for the trapped 
particle component (blue lines), streaming particle component (red 
lines), and their sum (green lines). Solid lines show spectra from 
FLASH and dashed lines represent spectra from AGILE-IDSA. Left 
panels show neutrino spectra and right panels show anti-neutrino 
spectra. Each row indicates spectra at a different radius. 

the model sll.O (527.0) when the shock reaches the Si/O 
interface. 

Figure 7 shows the time evolution of electron-type neu¬ 
trino luminosity and mean neutrino energy for the four 
considered models. The values are sampled at a radius 
of 500 km in both codes. The electron neutrino lumi¬ 
nosities and electron anti-neutrino luminosities are sim¬ 
ilar for the two codes around bounce and early post¬ 
bounce. However, at ^ 100 ms post bounce, the neu¬ 
trino and anti-neutrino luminosities in FLASH become 
slightly higher than in AGILE-IDSA for the progenitors 
515.0, 521.0 and 527.0. The maximum difference is an 
increase of about < 15% in FLASH at ^ 150 ms post 
bounce. After ~ 200 ms post bounce, the neutrino lu¬ 
minosities are back to the same value in both codes but 
the anti-neutrino luminosity remains slightly higher in 
FLASH simulations. For the progenitor 511.0, the elec¬ 
tron neutrino luminosity is the same in both codes but 
the electron antineutrino luminosity is slightly higher in 
FLASH. However, we note that the mean neutrino and 
anti-neutrino energies in FLASH are about 10 — 20% lower 
than in AGILE-IDSA. This difference could originate from 
the lower neutrino fraction in the low density region in 
FLASH (see Figure 3), since we measure the neutrino mean 
energy at r = 500 km. 

The kernel of our IDSA solver in FLASH is the same as 
the solver in AGILE-IDSA. The main differences are from 
the hydrodynamics, the explicit solver for the diffusion 
scalar a in Equation 10, the detailed treatment of the 
EOS, and potentially, also the gravity solver. The low- 
density extension for the EOS in FLASH and differences 
in the internal energy shift may also provide slightly dif¬ 
ferences. In principle, we should expect identical results 
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Figure 6. Time evolution of the shock radius for the four investi¬ 
gated progenitor models in both FLASH (solid lines) and AGILE-IDSA 
(dashed lines). Different colors represent different progenitor mod¬ 
els. All results are ID with LS220 EOS and without NES. The 
differences in the models sll.O and s27.0 after 0.1 s are due to a 
different handling of shell interfaces and a large diffusivity of the 
implicit adaptive mesh in AGILE-IDSA. 


Figure 5. Similar to Figure 2 but for the postbounce evolution. 

in both codes in ID. As presented above, although some 
little differences have been observed, most features are 
consistent and in nice agreement. Liebendorfer et al. 
(2005) have performed a code comparison of the Boltz¬ 
mann solver AGILE-BOLTZTRAN with the variable Edding¬ 
ton factor VERTEX code. Both codes have sophisticated 
physics input in spherical symmetry but different imple¬ 
mentations. As pointed out in Liebendorfer et al. (2005), 
the different grids in the Lagrangian or Eulerian coordi¬ 
nates produce late time differences in the shock evolu¬ 
tion when the shock runs through shell interfaces. Since 
FLASH is also an Eulerian code, similar differences be¬ 
tween FLASH and AGILE-IDSA as we have found here can 
be expected. 

3.2. Code Comparisons for the 15 M 0 and 20 M 0 
progenitors 

The si5 and s20 progenitors from Woosley & Heger 
(2007) have been widely studied in the literature, e.g. by 
Bruenn et al. (2013), Dolence et al. (2015), Suwa et al. 
(2014), Hanke (2014, PhD Thesis), and Melson et al. 
(2015). In addition, Lentz et al. (2012a, b) have shown 
that different neutrino opacities and approximations 
could lead to different post-bounce evolutions by us¬ 
ing spherically symmetric AGILE-BOLTZTRAN simulations 
with the si5 progenitor. Therefore these two progen¬ 
itor models are very suitable for comparisons among 
different SN codes. It is generally agreed that the 
Vertex-Prometheus (Rampp & Janka 2002) and the 
CHIMERA (Bruenn et al. 2013) codes use state-of-the-art 
physics for neutrino transport. Therefore, it is worth 
to perform simulations with these progenitor models to 
have a direct comparison. However, there is still different 
physics employed in the above mentioned works. For in¬ 
stance, CASTRO simulations (Dolence et al. 2015), ZEUS 
simulations (Suwa et al. 2014), and our FLASH-IDS A 
simulations are Newtonian, but AGILE-BOLTZTRAN sim- 




Figure 7. The electron neutrino (left) and electron antineutrino 
(right) luminosities and mean energies as a function of time. The 
quantities are averaged in a 5 ms time interval. Different colors 
represent the different progenitor models. Solid lines represent the 
FLASH simulations and dashed lines indicate the AGILE-IDSA simu¬ 
lations. 

ulations (Lentz et al. 2012a,b), CHIMERA simulations 
(Bruenn et al. 2013), and Vertex-Prometheus simula¬ 
tions (Hanke, 2014, and Melson et al. 2015) are GR or 
post-Newtonian. Furthermore, Dolence et al. (2015) use 
Shen’s EOS while other groups use the LS220 EOS. 
Therefore in this Section, we only give a qualitative dis¬ 
cussion. 

To evaluate that the PD scheme could effec¬ 
tively represent NES, we perform three Newtonian 






























Table 3 

ID Simulation Results During Core Collapse 


Property 




Models 





sll.O 

sl5.0 

s21.0 

s27.0 

Progenitor mass (Mq) 


11 


15 


21 


27 

Progenitor compactness (£ 1 . 75 ) 

0.066+ 

0.54 

0.68 

0.53 

Progenitor compactness (£2.5) 

0.004 

0.15 

0.22 

0.23 


w PD 

wo PD 

w PD 

wo P D 

w PD 

wo P D 

w PD 

wo P D 

Abbreviation (1D-) 

DP11 

DA11 

DP15 

DA15 

DP21 

DA21 

DP27 

DA27 

Time at p c = 10 n g cm -3 (ms) 

176.3 

194.1 

267.6 

278.2 

269.7 

284.3 

273.0 

287.2 

Time at p c = 10 12 g cm -3 (ms) 

197.6 

224.2 

285.7 

300.8 

286.3 

304.0 

289.4 

306.9 

Time at p c = 10 13 g cm -3 (ms) 

201.9 

229.3 

289.6 

305.5 

290.2 

308.4 

293.2 

311.4 

Time at p c = 10 14 g cm -3 (ms) 

203.0 

230.7 

290.7 

306.8 

291.2 

309.7 

294.3 

312.7 

Time at bounce (ms) 

203.5 

231.3 

291.2 

307.3 

291.7 

310.2 

294.7 

313.2 

Bounce, central p (10 14 g cm -3 ) 

3.21 

3.39 

3.14 

3.43 

3.08 

3.35 

3.10 

3.45 

Bounce, central Y e 

0.287 

0.316 

0.282 

0.312 

0.279 

0.310 

0.279 

0.311 

Bounce, shock position (Mq) 

0.54 

0.77 

0.56 

0.76 

0.56 

0.75 

0.56 

0.75 


t Note that the compactness parameter of the sll.O progenitor is about ten times smaller than others due to a smaller core mass and 
light envelope. 



Figure 8. Time evolution of shock radius (left) and electron-type neutrino luminosities (right). Different colors represent different models 
from FLASH and AGILE-BOLTZTRAN with different neutrino physics. See Section 3.2 for detailed description. 


AGILE-BOLTZTRAN simulations with the sl5 progenitor 
and the LS220 EOS. The first simulation (model AB- 
NR) includes the NES; the second simulation (model AB- 
NR-NoNES) ignores the NES; and the third simulation 
(model AB-NR-Mix) includes NES only in the collapse 
phase. Figure 8 shows the shock-radius and neutrino lu¬ 
minosity evolutions of these three simulations together 
with FLASH-IDS A simulations with and without PD. Ig¬ 
noring NES (models FLASH-IDSA and AB-NoNES in 
Figure 8) gives a short-period of shock expansion at ~ 
10 ms postbounce. A signature from this can also be seen 
in the electron anti-neutrino luminosity. This behavior is 
consistent with the model “N-ReduceOp” in Lentz et al. 
(2012a), and the shock and electron-type neutrino lumi¬ 
nosity evolutions of our FLASH-IDSA simulation without 
PD is also consistent with our AGILE-BOLTZTRAN simu¬ 
lation (AB-NR-NoNES). 

The model AB-NR-Mix is nearly identical to the model 


AB-NR, demonstrating that NES is important mainly 
during the collapse phase (Figure 8). Fischer et al. 
(2012) showed that that NES plays a role in the late 
time PNS cooling and deleptonization (> Is) after 
the explosion. However, our simulations focus only on 
the first few hundred milliseconds. The PD scheme in 
model FLASH-PD greatly improves the post-bounce sim¬ 
ulation, making our FLASH-IDSA simulations closer to 
the AGILE-BOLTZTRAN simulation with full neutrino re¬ 
actions. Although there is no perfect match, we con¬ 
clude that the use of the PD scheme in the collapse phase 
could effectively take into account the NES. It should be 
noted that heavy neutrinos are treated by a simple leak¬ 
age scheme which could lead to some difference in our 
IDSA simulations as well. 

In Figure 9, we show the shock-radius evolution and 
neutrino signatures of our 2D FLASH-IDSA simulations 
of the si5 and s20 progenitors from Woosley & Heger 
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Figure 9. Time evolution of shock radius (left) and electron-type neutrino signals (right). Different colors represent different progenitor 
models. In the left panel, thick lines show the evolution of 2D simulations and thin lines represent ID simulations. In the right panel, 
thick lines indicate electron-type neutrino luminosities and thin lines show the RMS neutrino mean energies. Solid lines represent electron 
neutrino signals and dashed lines represent electron antineutrino signals. 


(2007). To get a fair comparison, we enable the PD 
scheme during collapse and use the LS220 EOS. Both 2D 
models explode at ~ 300 ms postbounce (see Table 4). 
The explosion time, £400 is defined by the time when 
the averaged shock radius exceeds 400 km and never re¬ 
cedes at the end of the simulation. Overall, the shock- 
radius evolutions are similar to the models B15-WH07 
and B20-WH07 in Bruenn et al. (2013) except that the 
CHIMERA simulations show earlier explosions. It should 
be noted that we use the old neutrino interactions from 
Bruenn (1985) and that our simulations are Newtonian. 
Muller et al. (2012) have shown that GR effects could 
enlarge the neutrino luminosities and therefore make it 
easier to explode. The model s20-2007 in Hanke et 
al. (2014) and Melson et al. (2015) shows a similar ex¬ 
plosion time as our model s20, but the shock radius at 
^150 ms shrinks to ^ 150 km due to GR effects. The 
reasons for the differences between Vertex-Prometheus 
and CHIMERA are still unclear, but the overall features 
for the progenitor s20 are still rather similar. However, 
the model sl5-2007 in Hanke et al. (2014) shows a very 
different result. The shock stalls for ^ 500 ms and then 
explodes at around ^ 600 ms. 

On the other hand, 2D CASTRO and ZEUS Newto¬ 
nian simulations by Dolence et al. (2015) and Suwa et al. 
(2014) did not obtain explosion with the sl5 and s20 pro¬ 
genitors. Our 2D simulations show a fast shock expan¬ 
sion after the prompt convection (~ 20 ms, see Figure 9). 
This is similar to what was observed in Dolence et al. 
(2015) but somewhat less dramatic. The prompt convec¬ 
tion and fast shock expansion coincide with an oscillation 
of the electron antineutrino luminosity at 10 — 20 ms (see 
Figure 9 and Figure 6 of Dolence et al. 2015). These 
could be caused by the reduced opacity or incomplete 
neutrino interactions as discussed before. Note that 
Dolence et al. (2015) use the Shen EOS which is con¬ 
sidered more difficult to lead to explosions than LS220 
(Couch 2013; Suwa et al. 2013). 


Suwa et al. (2014) also use the IDSA (without PD) 
but with spherical coordinates and the “Ray-by-ray” ap¬ 
proach. In principle, we should expect similar results, 
but the non-explosion of sl5 and s20 in Suwa et al. (2014) 
suggest that the different hydrodynamics code, geom¬ 
etry, resolutions, and multidimensional neutrino trans¬ 
port approximation may also cause significant differ¬ 
ences. Suwa et al. (2014) use 300 logarithmically spaced 
radial zones (from 1km to 5,000 km) and 1.4° angular 
resolution. This is roughly three times lower than our 
simulations. A detailed code comparison is therefore nec¬ 
essary. 

4. MULTI-PROGENITOR STUDY 

We perform ID and 2D simulations with si 1.0, sl5.0, 
s21.0 and s27.0 progenitor models from Woosley et al. 
(2002). Simulations start from the prebounce core col¬ 
lapse to several hundred milliseconds postbounce with 
and without the PD in the collapse phase. The former 
is important to effectively take NES into account. Ta¬ 
ble 3 shows the core properties of these four progenitors 
during collapse based on ID simulations. A summary 
of all performed simulations is shown in Table 4. The 
model abbreviations in Tables 3 and 4 are defined by 
a set of letters and numbers: The first two characters 
define the dimension of the model (ID or 2D); the first 
letter after the hyphen denotes the EOS of the model (L 
for LS220 and D for DD2); the second letter shows the 
transport scheme during the collapse (A for IDSA and P 
for PD); and the last two numbers specify the mass of 
the investigated progenitor model. A “’-07” at the end 
shows progenitor models from Woosley & Heger (2007), 
otherwise they are from Woosley et al. (2002). For in¬ 
stance, model 1D-DA15 means a ID simulation of the 
sl5.0 progenitor with DD2 EOS and using the IDSA in 
the collapse phase (i.e. effectively without NES). When 
we refer to models DA, we consider all models with “DA” 
in its abbreviations. 
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Table 4 

Simulation Parameters and Results 


Progenitor 

Abbreviation^ 

Af'r a 

1 v 1 collapse 

EOS 

^bounce 

(ms) 

*400 C 
(ms) 

*end 

(ms) 

*^diag 

(B) 

*WpNS F 

(M G ) 

-RpNS G 

(km) 


ID 

sl5 (W07) 

1D-LA15-07 

IDSA 

LS220 

237 

— 

763 

— 

1.85 

29.1 

sl5 (W07) 

1D-LP15-07 

PD 

LS220 

249 

— 

523 

— 

1.77 

36.1 

s20 (W07) 

1D-LP20-07 

PD 

LS220 

322 

— 

678 

— 

1.99 

30.4 

sll.O (W02) 

1D-LA11 

IDSA 

LS220 

206 

— 

794 

— 

1.48 

27.5 

sl5.0 (W02) 

1D-LA15 

IDSA 

LS220 

273 

— 

727 

— 

1.84 

29.5 

s21.0 (W02) 

1D-LA21 

IDSA 

LS220 

274 

— 

726 

— 

1.98 

30.4 

s27.0 (W02) 

1D-LA27 

IDSA 

LS220 

283 

— 

717 

— 

1.81 

31.2 

sll.O (W02) 

1D-DA11 

IDSA 

DD2 

231 

— 

766 

— 

1.47 

32.2 

sl5.0 (W02) 

1D-DA15 

IDSA 

DD2 

307 

— 

693 

— 

1.84 

34.6 

s21.0 (W02) 

1D-DA21 

IDSA 

DD2 

310 

— 

604 

— 

1.95 

37.1 

s27.0 (W02) 

1D-DA27 

IDSA 

DD2 

313 

— 

687 

— 

1.81 

36.1 

sll.O (W02) 

1D-DP11 

PD 

DD2 

203 

— 

761 

— 

1.47 

33.6 

sl5.0 (W02) 

1D-DP15 

PD 

DD2 

291 

— 

709 

— 

1.84 

36.6 

s21.0 (W02) 

1D-DP21 

PD 

DD2 

292 

— 

708 

— 

1.98 

36.6 

s27.0 (W02) 

1D-DP27 

PD 

DD2 

295 

— 

705 

— 

1.81 

37.1 

2D 

sl5 (W07) 

2D-LP15-07 

PD 

LS220 

249 

312 

524 

0.298 

1.69 

37.1 

s20 (W07) 

2D-LP20-07 

PD 

LS220 

324 

284 

490 

0.347 

1.86 

40.5 

sl5.0 (W02) 

2D-LA15 

IDSA 

LS220 

274 

209 

311 

0.464 

1.60 

46.8 

sl5.0 (W02) 

2D-LA151owtt 

IDSA 

LS220 

274 

210 

369 

0.523 

1.62 

43.5 

sll.O (W02) 

2D-DA11 

IDSA 

DD2 

232 

86 

374 

0.821 

1.31 

42.3 

sl5.0 (W02) 

2D-DA15 

IDSA 

DD2 

308 

186 

417 

0.506 

1.65 

43.5 

s21.0 (W02) 

2D-DA21 

IDSA 

DD2 

311 

189 

411 

0.635 

1.75 

44.8 

s27.0 (W02) 

2D-DA27 

IDSA 

DD2 

314 

162 

429 

0.248 

1.66 

42.3 

sll.O (W02) 

2D-DP11 

PD 

DD2 

204 

170 

376 

0.267 

1.37 

47.4 

sl5.0 (W02) 

2D-DP15 

PD 

DD2 

292 

275 

456 

0.255 

1.71 

45.4 

s21.0 (W02) 

2D-DP21 

PD 

DD2 

292 

282 

484 

0.417 

1.82 

44.8 

s27.0 (W02) 

2D-DP27 

PD 

DD2 

295 

199 

484 

0.157 

1.70 

43.5 


A Neutrino transport scheme during the collapse phase 
E Bounce time 

^ Explosion time after bounce 
D Termination time after bounce 

E Diagnostic explosion energy at the end of the simulation 
F PNS mass at the end of the simulation 

G PNS radius (determined at the average radius corresponding to p = 10 11 g cm -3 ) at the end of the simulation 
i The model abbreviations are defined by the model’s dimensionality, EOS, neutrino transport scheme, and progenitor mass. See 
Section 4 for a detailed description. 

^ The effective angular resolution in this model is a factor of two lower. 


4.1. Stellar Collapse and Core Bounce 

Simulations are started from the non-rotating, 
solar-metallicity pre-supernova progenitors from 
Woosley et al. (2002) without artificial perturbation. 
Couch & Ott (2013) and Mueller & Janka (2014) show 
that small perturbations on the Si/O interface during 
collapse could amplify post-shock turbulence and turn 
a model that failed to explode towards a successful 
explosion. In addition, Couch et al. (2015) performed 
3D simulations of the final few minutes of iron core 
growth in a massive star that includes Si shell burning. 
The results suggest that the non-spherical progenitor 
structure may have a strong impact on neutrino-driven 
explosions. In our models, we do not include these 
non-spherical features, and therefore the 2D simula¬ 
tions during collapse show nice agreement with ID 
simulations in all models. However, we will show that 
spherical variations due to different levels of electron 
deleptonization or neutrino reactions during collapse 
may also have a significant impact on neutrino-driven 
explosions. 

Figures 10-13 show the radial density, electron fraction, 


entropy, and radial velocity evolutions of 2D models with 
HS(DD2) EOS based on the IDS A or the PD at differ¬ 
ent time during the collapse phase. The ID data are not 
shown because they are barely distinguishable from the 
2D data before bounce. In the bounce profiles, the mod¬ 
els DP show a slightly lower central density but a higher 
density distribution outside of the iron core, due to an 
earlier bounce time. 

The core evolutions of different progenitors behave 
qualitatively similar during collapse, but show quanti¬ 
tive differences for the central electron fraction, central 
density, and bounce shock locations. Table 3 summa¬ 
rizes important quantities for the four different progen¬ 
itors during the collapse phase. Due to a higher initial 
central density of the progenitor s 11.0, the central den¬ 
sity of s 11.0 reaches p c = 10 11 g cm -3 ^100 ms earlier 
than other progenitors. Once the core densities in dif¬ 
ferent progenitors reach the same p c , the more massive 
progenitors collapse faster than other progenitors. 

In addition, since the electrons are highly degener¬ 
ate and can only gain energy, this means that neutri¬ 
nos lose their energy through NES and therefore es- 
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Figure 10. Angle averaged radial profiles of progenitor si 1.0 at 
different times during core collapse. All models are using the DD2 
EOS. Different colors indicate different averaged profiles at certain 
central densities and bounce. The solid lines show simulations with 
the PD, and the dashed line show simulations with the IDSA. 

cape more easily, accelerating the collapse process (Bethe 
1990). Therefore models DP (with effective NES) col¬ 
lapse faster than models DA (without NES) and have a 
lower Y e (more e— captures), Yi (-\-v— escape) and p c at 
core bounce. 

Core bounce is defined here by the first time when the 
maximum density in the core exceeds 2 x 10 14 g cm -3 and 
the maximum peak entropy is above 3 ke baryon -1 . At 
bounce, the bounce shock emerges at ~ O.55M 0 (defined 
by the enclosed mass within the shock front at bounce) 
in models DP and at ^ O.75M 0 in models DA. Since the 
core mass is proportional to Y e 2 (Yahil & Lattimer 1982), 
models DA have higher core mass than models DP at 
bounce. The highest infall velocity at bounce is about 
^ 70,000 km s -1 in all progenitor models. It should be 
noted that the Si/O interface, which corresponds to the 
high entropy gradient region at about 100 — 400 km in 
Figures 10-13, is at different radii for models DA and 
DP due to different bounce times, since the bounce time 
directly sets the location of the Si/O interface. These 
difference will strongly influence the shock evolution after 
bounce. 


Figure 11. Similar to Figure 10 



Radius (km) 


but for progenitor si5.0. 



4.2. Shock Expansion and Instabilities 

Figure 14 shows the angle averaged shock radius as a 
function of time in our ID and 2D simulations with the 
DD2 EOS. The shock radius of a sample ray is defined 
here by the highest radius where the entropy is above 
<s m in = 4.5 ks baryon -1 (s m i n = 6 ks baryon -1 for 
progenitor models si 1.0 and s27.0) after the postbounce 
time was reached t p b = 50 ms. In 2D simulations, we 
sample 180 radial rays. Before 50 ms postbounce, the 
shock front is defined by the minimum radial infall ve¬ 
locity. 


Figure 12. Similar to Figure 10 but for progenitor s21.0. 

ID and 2D models behave very similarly in the first 
few milliseconds until the bounce shock passes through 
the neutrino sphere at ^ 10 ms. At that time, a 
prompt entropy- and electron-driven convection occurs 
and makes the 2D simulations different from ID. Later 
on, this prompt convection causes a fast shock expan¬ 
sion, and changes the shock radius from r s h ^100 km to 
r s h ~ 200 km at t p b ~ 20 ms (see Figure 14). While the 
prompt convection is caused by the negative entropy and 
Y e gradient, it should be noted that this early fast shock 
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Figure 13. Similar to Figure 10 but for progenitor s27.0. 

expansion during £ p b ^ 10 — 20 ms may be amplified by 
our incomplete neutrino interactions and the old neutrino 
opacities, and by the grid-effect in cylindrical coordinates 
in consideration with the multidimensional treatment of 
neutrino diffusion, since a Rayleigh-Taylor bubble along 
the diagonal axis is observed at that time. A similar 
fast shock expansion is observed in Dolence et al. (2015) 
as well, but is not obvious in any other 2D simulation 
with spherical coordinates. However, unlike the results 
in Dolence et al. (2015), the shock radius in our simu¬ 
lations shrinks back to r s h ^150 km after this prompt 
convection at £ p b ~ 50 ms. 

The prompt and late convection can be understood 
from the local stability analysis via the Ledoux criterion 
(Ledoux 1947), 



( 12 ) 


where we assume Y/ ^ Y e for simplicity. The Brunt- 
Vaisala (BV) frequency, oobv, describes the linear growth 
frequency for convection. Follow the definition of 
(Buras et al. 2006; Ott et al. 2013), one obtains 


lobv = sign (Cl) 



(13) 


where <f> is the local gravitational potential and the ap¬ 
proximation d^/dr rsj —GM(r)/r ~ 2 was used. Once con¬ 
vection is active, another useful quantity to describe the 
strength of convention is the anisotropic velocity. We fol¬ 
low the definition from Takiwaki et al. (2012) and define 
the anisotropic velocity as 


^aniso 



(v r ~ (Vr) 4 J 2 + v l 

^ 4-7T 

(p) 4tt 


(14) 


In Figure 15, we show the evolution of the angle- 
averaged BV frequency and anisotropic velocity of our 
models DP. The BV frequencies are positive and high 
when the shock break through the neutrinosphere and 
the prompt convection happens at ^ 20 ms. Starting 
from the core regions, the anisotropic velocities become 
very strong after the prompt convection has developed, 
and therefore drive the fast-shock expansion at ^ 20 ms. 
However, the convection stops and the anisotropic ve¬ 
locity returns to small values at ~ 50 ms. Models DA 
behave similar but on a different time scale due to dif¬ 
ferent shock evolutions. We note that our estimate of 
Cl may be incorrect at small radii, where neutrinos are 
trapped, due to the ignorance of the contribution from 
neutrinos in Yj, but the overall features should be the 
same. 

For a given progenitor model, the bounce shocks have 
a similar shock strength in both models DA and DP. The 
difference of shock radius in models DA and DP during 
the prompt convection is mainly due to the difference of 
the bounce time, since different bounce time leads to a 
different shock strength when the bounce shock reaches 
the location of the negative entropy gradient. Since con¬ 
vection is suppressed in ID, we don’t observe this prompt 
convection in ID models. However, the prompt convec¬ 
tion in 2D models leads to a larger shock radius than in 
ID models at £ p b ^ 50 ms. 

After t p b ~ 50 ms, the shock stalls at r s h ^ 150 — 
200 km for another ^ 50 — 100 ms (except for model 2D- 
DAll). It is no surprise that all ID models fail to ex¬ 
plode due to the lack of convection. However, at first it 
seems surprising that all our 2D models explode within ^ 
100 — 300 milliseconds postbounce regardless of the pro¬ 
genitor mass and the inclusion of NES during collapse. 
But on second thought, one has to consider that our sim¬ 
ulations are Newtonian (comparing with Bruenn et al. 
2001; Liebendorfer et al. 2001; Muller et al. 2012) and 
do not include the most recent electron capture rates 
(Langanke et al. 2003; Hix et al. 2003). In general, New¬ 
tonian simulations disfavor explosions compared to gen¬ 
eral relativistic simulations (Muller et al. 2012), and also 
the old neutrino interactions and opacities could make 
our simulations too optimistic with respect to explosions. 
Furthermore, we are using the DD2 EOS, which has more 
realistic nuclear matter properties than the LS220 EOS 
that was used in the above mentioned simulations. In 
Appendix B, we show that DD2 leads in our simula¬ 
tions to more favorable conditions for explosions than 
LS220. Table 4 summarizes the essential parameters of 
our ID and 2D simulations at the end of the simulations. 
Model 2D-DA11 is the fastest explosion model that ex¬ 
plodes already at £400 = 86 ms by neutrino-driven con¬ 
vection. For the other DA models, the explosion time 
is £400 ^ 160 — 190 ms. In addition, because of the ef¬ 
fective NES, models DP evolve more slowly and explode 
at a late £400 ^ 200 — 280 ms (£400 = 170 ms for model 
2D-DP11, see Table 4 for detailed information). 

In model 2D-DP11, the shock radius first expands 
to r s h > 300 km at ~ 100 ms postbounce when it 
reaches the Si/O interface, but temporarily drops back to 
v s h < 300 km at ~ 130 ms postbounce before it explodes. 
Progenitor s 15.0 and s21.0 explode later than progeni¬ 
tor s 11.0 and s27.0 but at a similar time in both models 
2D-DA and 2D-DP. In Figures 16-19, we compare the 
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Figure 14. Average shock radius vs. postbounce time for different progenitors and neutrino transport approximations. All models are 
using the DD2 EOS. Different colors indicate different progenitor models. The left panel represents simulations without the PD, and the 
right panel represents simulation with an effective inclusion of NES by the PD. The solid lines show 2D simulations, and the dashed lines 
show ID simulations. 


average radial profiles of density, entropy, electron frac¬ 
tion, and radial velocity at 1, 3, 50, 100, and 200 ms 
postbounce for models DA and DP. We terminate the 
simulations at ^ 700 — 800 ms postbounce in ID, and 
~ 300 — 500 ms postbounce in 2D. 

The 2D entropy distribution of models 2D-DP are 
shown in Figure 20 at 150, 250, and 300 ms postbounce. 
At ~ 150 ms, the convection is getting stronger (see Fig¬ 
ure 15) but the shock radius distribution is still spheri¬ 
cally symmetric. Later on, at ^ 200 ms the SASI starts 
in models 2D-DP15, 21, and 27 but the models of pro¬ 
genitor sll.O show mainly convection without SASI. The 
SASI activities can be seen in Figure 21, where we show 
the normalized coefficients ai of the decomposition of 
the shock radius r s h($) into Legendre polynomials Pi 
(Muller et al. 2012). ai can be calculated by: 

2/ + l r 

ai = —-— J r sh (0)Pid cos 0, (15) 

and ao corresponds to the averaged shock radius. For the 
progenitor si 1.0, there is no obvious evidence of SASI 
activities in both models 2D-DA11 and 2D-DP11. The 
amplitude of the normalized coefficients for l = 1, 2 and 3 
modes are small and within the same order of magnitude. 
For progenitors sl5.0, s21.0, and s27.0, the SASI activ¬ 
ities could be seen and start to grow at ~ 200 ms post¬ 
bounce. After ~ 200 ms postbounce, the dipole (1 = 1) 
and quadrupole (l = 2) modes grow to ai/ao ^ 0.2 in 
the progenitors si5.0 and s21.0, and ai/ao ~ 0.1 in pro¬ 
genitor s27.0. The amplitudes do not show significant 
differences between models DA and DP, but the starting 
time of the high growth rates of the amplitudes corre¬ 
sponds to the time of fast shock expansion. SASI ac¬ 
tivities could also be seen in Figure 22 for the entropy 
distribution along the north and south poles. 

To verify the code convergence, we have also performed 
a low-resolution run by reducing the angular resolution 
by a factor of 2 (model 2D-LA151ow). We find that the 


explosion time £400 is 1ms delayed and the shock expan¬ 
sion evolves slightly slower than for the standard resolu¬ 
tion run (model 2D-LA15). Overall, we find no signifi¬ 
cant differences between the low-resolution run and the 
standard run. 

4.3. Neutrino Heating and Explosion Energy 

The energy released in neutrinos is about 10 53 erg in 
a typical CCSN (Colgate & White 1966). To power the 
observed kinetic energy of a CCSN (^ 10 51 erg = 15), 
the baryonic matter has to absorb ~ 1% of the neutrino 
energy (Bethe & Pizzochero 1990). However, the explo¬ 
sion energies in most published 2D models are still lower 
than the standard IB (except for those of Bruenn et al. 
2013, 2014). The real explosion energy is not straight¬ 
forward to calculate and the definition may differ from 
group to group. Figure 23 shows the “diagnostic energy” 
of our models 2D-DA and 2D-DP. The diagnostic energy, 
5dia, is defined by 

Edm = f pf-HotdV, (16) 

d&tot >0 

where the volume integration is performed over the re¬ 
gion where the total specific energy, e to t, 

e to t = (e- e 0 ) + h 2 + $ (17) 

is positive, e represents the specific internal energy (ther¬ 
mal energy plus binding energy), and \v 2 and <f> are the 
specific kinetic and gravitational energies, respectively, 
eo is the reference energy value that is defined by the 
minimum specific internal energy at the beginning of the 
simulation, which leads to a negligible diagnostic explo¬ 
sion energy at beginning. We have checked that alterna¬ 
tive values for the reference energy, that are suggested 
in the literature, do not have a significant effect on our 
final results. Depending on progenitor and simulation 
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Figure 15. The color maps show the time evolution of the angle-averaged Brunt-Vaisala (BV) frequency ujbv in units of ms -1 (left) and 
the anisotropic velocity v an i so in units of the speed of light c (right) in our models 2D-DP. The red lines show the time evolution of the 
maximum shock radius and the black (left) and white (right) lines show the time evolution of the averaged shock radius. 


domain, eo is about ~ 10 17 erg g _1 . We find that the 
diagnostic energy is insensitive to eo if we make it a few 
times larger or smaller. Note that this expression does 
not take into account properly the different binding en¬ 
ergy contributions from different nuclei to thermal en¬ 
ergy. Therefore the so-called recombination energy is not 
included, e corresponds only to the thermal energy that 
actually should be used in 5di a , if the composition was 
dominated by 12 ( 7 , or if the composition has a similar 
average binding energy like 12 C. The diagnostic explo¬ 


sion energies (Figure 23) of our models 2D-DP are be¬ 
tween 0.1 B and 0.45 at ^ 400 ms postbounce and are 
still increasing at the end of our simulations. The mod¬ 
els 2D-DA, which ignored the NES in the collapse phase, 
have a much stronger explosion energy, 5 exp ~ 0.2 — 0.85 
at ~ 400 ms postbounce. This result is consistent with 
the ID results in Bruenn (1989). It should be noted 
that the only difference between models 2D-DA and 2D- 
DP is the deleptonization method during collapse. The 
postbounce physics employed are identical in both mod- 
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Figure 16. Angle averaged radial profiles of progenitor si 1.0 at 
different postbounce times. All models are evolved in 2D and are 
using the DD2 EOS. Different colors indicate different postbounce 
times. The solid lines show simulations with the PD, and the 
dashed lines show simulations without the PD. 


ml5 



Figure 17. Similar to Figure 16 but for progenitor sl5.0. 

els 2D-DA and 2D-DP, suggesting that the initial condi¬ 
tions at bounce may dramatically affect the postbounce 
evolutions. Furthermore, note that it is difficult to com¬ 
pare the results of models DP with other groups because 
there are differences in the input physics and methods. 
Additionally, we show the diagnostic explosion energies 
of models 2D-LA15 and 2D-LA151ow for a comparison 



Figure 19. Similar to Figure 16 but for progenitor s27.0. 

with LS220 EOS and with low resolution. We find that 
the model with LS220 EOS leads to a slight lower diag¬ 
nostic explosion energy due to a delay of the explosion. 
The low resolution model also gives a similar results, sug¬ 
gesting that the fast explosion in our models is not due 
to insufficient resolution (Figure 23). 

Ugliano et al. (2012); Nakamura et al. (2014), and 
Perego et al. (2015) suggest that there is a correlation 
of the compactness parameter with the explosion en¬ 
ergy. We do see this trend in our simulations ex- 
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Figure 20. The entropy distributions of our 2D models with the 
PD and the DD2 EOS (See Table 4). Each frame shows a section of 
the domain spanning 800 km. The color scale indicates the entropy 
in ^B/baryon. 

cept for the progenitor s 11.0, if we define the com¬ 
pactness parameter at an enclosed mass of 1.75M 0 , i.e. 

6-75 = fl(M=i.75Me)/i000km (O’Connor & Ott 2011). In 
our simulations, the progenitor s 11.0, which has the low¬ 
est compactness parameter, has the second highest diag¬ 
nostic explosion energy (see Tables 3 and 4). However 
this trend in our simulations will disappear if we use 
<^ 2 . 5 - It should be noted that we only have four progeni¬ 
tor models in our simulations, and the correlation found 
in Ugliano et al. (2012) and Perego et al. (2015) show 
a huge scatter between different models, indicating that 






Time (sec) 
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Figure 21. First, second, and third coefficients ai, a 2 , and <23 of 
the spherical decomposition of the shock radius into Legendre poly¬ 
nomials, normalized to the average shock radius (no) for different 
progenitors. 


the compactness parameter oversimplifies the complexity 
of the explosion mechanism. 

Figure 24 shows the neutrino luminosity and mean 
energy of electron type neutrinos and anti-neutrinos as 
functions of time. The neutrino mean energy, E v is de¬ 
fined by E v = where / is the distribution 

function. The neutrino luminosities in models DA after 
^ 50 ms grow faster than in models DP but the peak 
luminosities after ^100 ms are similar in both DA and 
DP models. The neutrino mean energies show a similar 
trend but the differences are less pronounced. We note 
that the neutrino luminosities and mean energies in our 
simulations are of a similar order of magnitude as the 
typical results in the literature (Janka 2012). Therefore, 
the high explosion energies in our simulations require a 
larger mass in the gain region or a higher heating effi¬ 
ciency than other investigators. 

Figure 25 presents the 2D net heating and cooling 
rates in models 2D-DP. The red color shows the heat¬ 
ing region and the blue color indicates the cooling re¬ 
gion. The PNS radius (defined by the average radius 
of density p = 10 11 g cm -3 ) is plotted as well in Fig¬ 
ure 25. The gain region is defined as the post-shock re¬ 
gion with positive net heating. In Figure 26, we show 
the integrated quantities of net heating (Q n et)> total 
mass in the gain region, and the heating efficiency as 
functions of time. The heating efficiency is defined by 
?7heat = ^5 net/(-^ zas, gain T A^ e g a i n ), since the heating orig¬ 
inates mainly from the electron type neutrinos. 

The early prompt convection shows a very high heat¬ 
ing efficiency peak at ~ 20 ms, which could be associated 
with a grid-effect in cylindrical coordinates. However, 
this prompt heating does not affect the explosion energy 
directly since matter is still bound at that time. The 
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Figure 22. Entropy along the north and the south poles as a function of time for different progenitor masses and for models DA (left) 
and DP (right) in Table 4. 


shock expansion and the early convection at t p b ^ 50 ms 
enlarge the gain radius and the mass in the gain region. 
After t p b ~ 100 — 200 ms, the mass in the gain region con¬ 
tinues to increase in time. We noted that the masses in 
the gain region O.IMq) in our simulations are higher 
than the values that are reported by other groups, ex¬ 
plaining the high diagnostic explosion energy in our sim¬ 
ulations. 

The models 2D-DA11 and 2D-DP11 have the lowest 
net heating and heating efficiency among the models (ex¬ 
cept during the prompt convection period), but have the 
second highest mass in the gain region and diagnostic 
energy. This is probably due to the fast shock expansion 
(see Figure 14) and low mass-accretion in the progenitor 
511.0. 

4.4. Proto Neutron Stars 

During its evolution, the PNS shrinks due to the neu¬ 
trino cooling behind the gain radius. In Figure 27, we 


show the mass and radius evolutions of the PNSs after 
bounce. An initial oscillation of the PNS is taking place 
at ^ 10 ms, when the shock breaks through the neutrino 
sphere. Later on, mass continuously accretes onto the 
PNS for several hundred milliseconds so that it reaches 
^ 1.3 — 2.0 Mq at the end of the simulations. At that 
time, the PNS radii have shrunken to ^ 40 km. The 
mass and radii of the PNSs at the end of the simulations 
are summarized in Table 4. 

The convection and SASI activities help the PNS to 
accrete mass during the setup of an explosion. Accre¬ 
tion funnels are created along the high entropy convec¬ 
tive plumes. This can be seen in Figure 20. Therefore, 
the progenitor s 11.0 which has the less SASI activities 
and accretion funnels, shows the lowest PNS mass incre¬ 
ment and the growth rate of the PNS mass is the lowest 
as well. Furthermore, the progenitor s 11.0 has the low¬ 
est density distribution outside the iron core, which gives 
the lowest mass accretion rate. After ^ 400 ms, the mass 
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Figure 23. The diagnostic explosion energies defined in Section 4 
as a function of time. Different colors represent the different pro¬ 
genitor models. Solid lines represent the models 2D-DP and dashed 
lines indicate the models 2D-DA in Table 4. The magenta and yel¬ 
low lines use LS220 EOS, while other lines use DD2 EOS. 
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Figure 24. Similar to Figure 7 but for the models DP (solid lines) 
and the models DA (dashed lines) in Table 4. 

increment of the PNS becomes small but the mass of the 
PNS is still increasing at the end of the simulations. We 
find that although the growth rate of the PNS mass in 
models 2D-DA and 2D-DP are different, the final PNS 
mass in models 2D-DP are only slightly larger (< 5%), 
and the radii of PNSs in models 2D-DA are slightly lower 
than that in models 2D-DP, giving a more compact PNS 
core in models 2D-DA. This can also be related to the 
earlier explosion times of models DA. 

Figure 28 shows the central density, central electron, 
lepton fraction, central entropy and central temperature 
as functions of time. Models 2D-DA show higher cen¬ 
tral densities and higher central Y e (higher electron pres¬ 
sures), but lower central entropies and central temper- 


Figure 25. The net heating (red) and cooling (blue) distributions 
of our models 2D-DP (See Table 4). The solid black lines show the 
radius of the PNS. 

atures than the models 2D-DP, giving a more compact 
PNS core. However, the density in the outer layer of 
the PNS is higher and more convective in models 2D-DP 
than that in models 2D-DA. Therefore, models 2D-DP 
have a higher PNS mass and radius (see Figure 27). 

5. CONCLUSIONS 

We have performed self-consistent CCSN simulations 
for the four nonrotating, solar-metallicity progenitors, 
sll.O, sl5.0, s21.0, and s27.0 from Woosley et al. (2002) 
and sl5 and s20 from Woosley & Heger (2007) using the 
AMR code FLASH with an IDS A for neutrino transport. 
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Figure 27. The time evolution of the PNS mass (left) and radius 
(right). The PNS mass and radius are determined at the average 
radius corresponding to p = 10 11 g cm -3 . Different colors repre¬ 
sent different progenitor stars. Solid lines show the models 2D-DP 
and dashed lines show the models 2D-DA in Table 4. 
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Figure 26. The integral quantities of net heating (top), mass 
in the gain region (middle), and heating efficiency (bottom) as 
functions of time. See Section 4 for the definition of these integral 
quantities. Different colors represent different progenitor models. 
Solid lines represent the models 2D-DP and dashed lines indicate 
the models 2D-DA in Table 4. 

A very good agreement of FLASH-IDSA with AGILE-IDSA 
is shown in Section 3.1, though some small differences 
still exist. In addition, we have provided a comparison of 
our multi-dimensional IDSA results for the sl5 and s20 
progenitors from Woosley & Heger (2007) with results 
in the existing literature. However, a detailed compari¬ 
son with an exchange of data and collaboration among 
groups would be necessary to understand remaining dif¬ 
ferent results with regard to the underlying physics. 
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Figure 28. The time evolution of central density (p c , upper left), 
central electron fraction ( Y e ) and lepton fraction (Yj, upper right), 
central entropy (lower left), and central temperature (lower right). 
Different colors represent different progenitor stars. Solid lines 
show the models 2D-DP and dashed lines show the models 20- 
DA in Table 4. In the upper right panel, the thin lines represent 
the central lepton fraction and the thick lines represent the central 
electron fraction. 

The new SN EOS HS (DD2) is employed in this study 
and a comparison with the standard LS220 EOS is dis¬ 
cussed in Appendix A. It is found that the DD2 EOS 
shows a faster shock expansion, and a higher mass in 
the gain region, while the neutrino luminosities are sim¬ 
ilar, causing the DD2 EOS to lead to slightly earlier and 
stronger explosions than the LS220. 

We have presented two sets of simulations which com¬ 
pared the neutrino transport with the IDSA (but ignor¬ 
ing NES; abbreviated DA) with the parametrized delep- 
tonization scheme (which includes NES effectively; ab¬ 
breviated DP) during collapse. The results show clearly 
that the treatment of neutrino weak interactions and the 
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level of deleptonziation during the collapse phase have a 
significant impact on the neutrino-driven explosions. All 
our 2D models explode within ~ 100 — 300 ms. Models 
without NES explode much easier, stronger and faster 
than models with NES, The diagnostic explosion energy 
in our simulations are around 0.1 — 0.5 B in models 2D- 
DP (see Table 4) and around 0.2 — 0.8 B in models 20- 
DA. Our explosion energies are likely to decrease when 
we include general relativistic effects and better electron- 
capture rates in future models. 

In our simulations, we do not use the typical RbR ap¬ 
proach for the neutrino transport. Instead, we solve 
the diffusion part in multi-dimensions, improving the 
neutrino transport in angular and temporal directions. 
With this approach, the prompt convection after the core 
bounce causes a fast shock expansion (t p b ~ 50 ms), en¬ 
larging the gain region at late time, and therefore helps 
to increase the explosion energy at late time. However, 
this prompt convection may be amplified by grid-effects 
in cylindrical coordinates. To address this issue, a full 
3D simulation will be necessary or a detailed comparison 
with simulations in spherical coordinates. Furthermore, 
we have found SASI activities at late time in progenitors 
515.0, s21.0, and 527.0, but in most cases, the explosions 
are mainly due to neutrino-driven convection. Future 
work will relax the constraints imposed by axisymme- 
try. Models in 3D will permit to study a more realistic 
turbulence cascade and rotation. 
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Figure 29. Similar to Figure 10 but for models 2D-DA15, and 
2D-LA15 in Table 4. 
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APPENDIX A 

A COMPARSION BETWEEN DD2 AND LS220 

Several simulations in this paper use the new DD2 
EOS, while most multi-dimensional simulations found in 
the literatures use LS220. We therefore briefly describe 
the qualitative differences between the simulations with 
the DD2 and the LS220 EOS in this section. We perform 
ID and 2D simulations of the progenitor 515.0 with DD2 
and LS220 EOS (models 1D-DA15, 2D-DA15, 1D-LA15 
and 2D-LA15 in Table 4). Hereafter, we use DA15 to 
refer to both models 1D-DA15 and 2D-DA15, and simi¬ 
larly we use LA15 to refer to the models 1D-LA15 and 


Figure 30. Similar to Figure 24 but for models 1D-DA15, 1D- 
LA15, 2D-DA15, and 2D-LA15 in Table 4. 

2D-LA15. Their main features are summarized in Ta¬ 
ble 4). The PD scheme is turned off in order to study 
the impact of the EOS on the neutrino transport. 

The angle-averaged radial profiles of density, electron 
fraction, entropy, and radial velocity are shown in Fig¬ 
ure 29 for the prebounce phase. As discussed in Sec¬ 
tion 4, the radial profiles in the collapse phase are nearly 
identical in both ID and 2D. In addition, the density and 
radial velocity evolutions in models DAI5 and LAI5 are 
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Figure 31. Similar to Figure 28 but for models 1D-DA15, 1D- 
LA15, 2D-DA15, and 2D-LA15 in Table 4. 
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Figure 32. Similar to Figure 27 but for models 2D-DA15, and 
2D-LA15 in Table 4. 

also similar during the collapse but models DAI5 collapse 
slower than models LAI5. Therefore, models DAI5 reach 
core bounce about ^ 30 ms later than models LAI5. 

While there are no significant differences of the neu¬ 
trino luminosity and mean energy between models DAI5 
and LA15 (Figure 30), models DA15 show a lower central 
density and higher PNS radius than LAI5. The time evo¬ 
lutions of central density and electron fraction of models 
DAI5 and LAI5 are shown in Figure 31. Furthermore, 
models DAI5 have a higher central electron and lepton 
fraction than models DL15. 2D models also show a lower 
PNS mass and a higher PNS radius than ID (Figure 32), 
but the PNS masses are similar in both DD2 and LS220 
(see Table 4). 

Since the bounce time in models DAI5 take longer than 
models LA15, the Si/O interface reaches a smaller radius 
in models DA15 (^ 90 km) than in LA15 (~ 120 km, 
see Figure 29). The bounce shock hits the Si/O inter¬ 
face earlier in models DAI5 than in LAI5, generating a 
prompt convection and fast shock expansion at ^ 20 ms 
postbounce. The earlier interaction between the bounce 
shock with the Si/O interface in model 2D-DA15 makes 
the shock radius in model 2D-DA15 larger than in sim¬ 
ulation 2D-LA15. Furthermore, the larger shock radius 
in model 2D-DA15 produces a higher mass in the gain 
region (see Figure 33), suggesting that the models with 
the DD2 EOS explodes more easily than the models with 
the LS220 EOS. The averaged shock radius evolution of 
models DA15 and LA15 can be seen in Figure 34. The 
explosion time (£ 400 ) in 2D-DA15 is about 23 ms earlier 
than in 2D-LA15. In Figure 35, we compare the radial 
density, electron fraction, entropy, and velocity profiles 


Figure 33. Similar to Figure 26 but for models 2D-DA15, and 
2D-LA15 in Table 4. 
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Figure 34. Similar to Figure 14 but for models 1D-DA15, 1D- 
LA15, 2D-DA15, and 2D-LA15 in Table 4. 


in ID and 2D for DD2 and LS220 at 1, 3, 50, 100, and 
200 ms postbounce. Before ^ 20 ms postbounce, models 
LSI5 show a faster shock expansion than models DAI5 
but the shock radius in DAI5 become larger after the 
Si/O interface has reached the shock at ^ 20 ms. 

Figure 36 shows the normalized decomposition into 
Legendre polynomials (Equation 15) of the shock radius 
of models 2D-DA15 and 2D-LA15. The corresponding 
entropy distribution at the north and south pole is shown 
in Figure 37. The SASI activities start at ^ 0.2 s in both 
2D-DA15 and 2D-LA15 after the Si/O interface reach the 
shock and the amplitude of the normalized coefficients 
ai/ao do not show significant difference. We remark that 
a more through investigation of the EOS effects will be 
the subject of a future study. 

REFERENCES 


Abdikamalov, E., Ott, C. D., Radice, D., et al. 2014, 
arXiv: 1409.7078 

Bethe, H. A. 1990, Reviews of Modern Physics, 62, 801 
Bethe, H. A., & Pizzochero, P. 1990, ApJ, 350, L33 
Berninger, H., Frenod, E., Gander, M., Liebendorfer, M., & 

Michaud, J. 2013, SIAM Journal on Mathematical Analysis, 45, 
3229 

















22 



Radius (km) Radius (km) 


Figure 35. Similar to Figure 16 but for models 1D-DA15, 1D- 
LA15, 2D-DA15, and 2D-LA15 in Table 4. 




Figure 36. Similar to Figure 21 but for models 2D-LA15 (Left) 
and 2D-LD15 (Right) in Table 4. 


Blondin, J. M., Mezzacappa, A., & DeMarino, C. 2003, ApJ, 584, 
971 

Brandt, T. D., Burrows, A., Ott, C. D., & Livne, E. 2011, ApJ, 
728, 8 

Bruenn, S. W. 1985, ApJS, 58, 771 

Bruenn, S. W. 1989, ApJ, 340, 955 

Bruenn, S. W. 1989, ApJ, 341, 385 

Bruenn, S. W., De Nisco, K. R., & Mezzacappa, A. 2001, ApJ, 
560, 326 

Bruenn, S. W., Mezzacappa, A., Hix, W. R., et al. 2013, ApJ, 
767, LL6 

Bruenn, S. W., Lentz, E. J., Hix, W. R., et al. 2014, 
arXiv: 1409.5779 

Buras, R., Rampp, M., Janka, H.-T., & Kifonidis, K. 2006, A&A, 
447, 1049 

Burrows, A. 2013, Reviews of Modern Physics, 85, 245 

Burrows, A., Young, T., Pinto, P., Eastman, R., & Thompson, 

T. A. 2000, ApJ, 539, 865 

Colella, P., & Woodward, P. R. 1984, Journal of Computational 
Physics, 54, 174 

Colgate, S. A., & White, R. H. 1966, ApJ, 143, 626 

Couch, S. M. 2013, ApJ, 765, 29 

Couch, S. M., Graziani, C., & Flocke, N. 2013, ApJ, 778, 181 

Couch, S. M., & Ott, C. D. 2013, ApJ, 778, LL7 


Couch, S. M., & O’Connor, E. P. 2014, ApJ, 785, 123 
Couch, S. M., Chatzopoulos, E., Arnett, W. D., & Timmes, F. X. 
2015, arXiv: 1503.02199 

Dolence, J. C., Burrows, A., & Zhang, W. 2015, ApJ, 800, 10 
Dubey, A., Reid, L. B., & Fisher, R. 2008, Physica Scripta 
Volume T, 132, 014046 

Fischer, T., Martfnez-Pinedo, G., Hempel, M., & Liebendorfer, 

M. 2012, Phys. Rev. D, 85, 083003 
Fischer, T., Hempel, M., Sagert, I., Suwa, Y., & Schaffner-Bielich, 
J. 2014, European Physical Journal A, 50, 46 
Foglizzo, T., Kazeroni, R., Guilet, J., et al. 2015, 
arXiv: 1501.01334 

Fryer, C., Benz, W., Herant, M., & Colgate, S. A. 1999, ApJ, 516, 
892 

Fryxell, B., Olson, K., Ricker, P., et al. 2000, ApJS, 131, 273 
Hanke, F., Marek, A., Muller, B., & Janka, H.-T. 2012, ApJ, 755, 
138 

Hannestad, S., & Raffelt, G. 1998, ApJ, 507, 339 
Hempel, M., & Schaffner-Bielich, J. 2010, Nuclear Physics A, 837, 
210 

Hempel, M., Fischer, T., Schaffner-Bielich, J., & Liebendorfer, M. 

2012, ApJ, 748, 70 
Hempel, M. 2014, arXiv: 1410.6337 

Hempel, M., Hagel, K., Natowitz, J., Ropke, G., & Typel, S. 

2015, arXiv: 1503.00518 

Hix, W. R., Messer, O. E., Mezzacappa, A., et al. 2003, Physical 
Review Letters, 91, 201102 

Janka, H.-T. 2012, Annual Review of Nuclear and Particle 
Science, 62, 407 

Janka, H.-T., & Mueller, E. 1996, A&A, 306, 167 
Kruger, T., Tews, I., Hebeler, K., & Schwenk, A. 2013, 

Phys. Rev. C, 88, 025802 & Liebendorfer, M. 2011, ApJS, 195, 
20 

Kuroda, T., Kotake, K., & Takiwaki, T. 2012, ApJ, 755, 11 
Kuroda, T., Takiwaki, T., & Kotake, K. 2015, arXiv: 1501.06330 
Langanke, K., Martfnez-Pinedo, G., Sampaio, J. M., et al. 2003, 
Physical Review Letters, 90, 241102 
Lattimer, J. M., Sz Douglas Swesty, F. 1991, Nuclear Physics A, 
535, 331 

Ledoux, P. 1947, ApJ, 105, 305 

Lee, D. 2013, Journal of Computational Physics, 243, 269 
Lentz, E. J., Mezzacappa, A., Messer, O. E. B., et al. 2012, ApJ, 
747, 73 

Lentz, E. J., Mezzacappa, A., Messer, O. E. B., Hix, W. R., & 
Bruenn, S. W. 2012, ApJ, 760, 94 
Liebendorfer, M., Mezzacappa, A., Sz Thielemann, F.-K. 2001, 
Phys. Rev. D, 63, 104003 

Liebendorfer, M., Mezzacappa, A., Thielemann, F.-K., et al. 

2001, Phys. Rev. D, 63, 103004 

Liebendorfer, M., Messer, O. E. B., Mezzacappa, A., et al. 2004, 
ApJS, 150, 263 

Liebendorfer, M. 2005, ApJ, 633, 1042 

Liebendorfer, M., Rampp, M., Janka, H.-T., Sz Mezzacappa, A. 
2005, ApJ, 620, 840 

Liebendorfer, M., Whitehouse, S. C., & Fischer, T. 2009, ApJ, 
698, 1174 

Livne, E., Burrows, A., Walder, R., Lichtenstadt, I., & 

Thompson, T. A. 2004, ApJ, 609, 277 
Melson, T., Janka, H.-T., & Marek, A. 2015, ApJ, 801, L24 
Melson, T., Janka, H.-T., Bollig, R., et al. 2015, ApJ, 808, L42 
Mezzacappa, A., Sz Bruenn, S. W. 1993, ApJ, 405, 669 
Mezzacappa, A., & Bruenn, S. W. 1993, ApJ, 410, 740 
Mueller, B., & Janka, H.-T. 2014, arXiv: 1409.4783 
Muller, B., Janka, H.-T., & Marek, A. 2012, ApJ, 756, 84 
Muller, B., Janka, H.-T., & Heger, A. 2012, ApJ, 761, 72 
Murphy, J. W., & Burrows, A. 2008, ApJ, 688, 1159 
Nakamura, K., Takiwaki, T., Kuroda, T., & Kotake, K. 2014, 
arXiv: 1406.2415 

Nomoto, K., & Hashimoto, M. 1988, Phys. Rep., 163, 13 
Obergaulinger, M., Janka, H.-T., & Aloy, M. A. 2014, MNRAS, 
445, 3169 

O’Connor, E., & Ott, C. D. 2010, Classical and Quantum 
Gravity, 27, 114103 

O’Connor, E., & Ott, C. D. 2011, ApJ, 730, 70 
O’Connor, E., & Ott, C. D. 2013, ApJ, 762, 126 
O’Connor, E. 2014, arXiv: 1411.7058 























23 



Time (ms) 


Time (ms) 


Figure 37. Similar to Figure 22 but for models 2D-LA15 (Left) and 2D-DA15 (Right) in Table 4. 
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